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We investigate the behavior of the periodic Anderson model in the presence of d-f Coulomb 
interaction (Udf). Three different metliods are applied: mean-field theory, variational calculation, 
and exact diagonalization of finite chains. The variational approach is based on the Gutzwiller trial 
wave function and we numerically optimize the variational parameters. Two quantum critical points 
(QCPs) exist for a critical value of Udf, where the valence susceptibility diverges. We derive the 
critical exponent for the valence susceptibility and investigate how the position of the QCP depends 
on the parameters of the Hamiltonian. For larger values of Udf, the Kondo regime is bounded by 
two first-order transitions. These first-order transitions merge into a triple point at a certain value 
of Udf ■ For even larger Udf valence skipping occurs. The results of the three methods are compared 
in detail. 

PACS numbers: 71.10.Fd, 71.27.-|-a, 75.30.Mb 
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I. INTRODUCTION 

Heavy-fermion compounds often show remarkable phe- 
nomena like unconventional superconductivity or un- 
usual Fermi-liquid state. It turned out that in those 
compounds, for instance Celnsji whose superconduct- 
ing state is unconventional, the pairing between electrons 
is mediated by antiferromagnetic spin fluctuations. The 
superconducting phase is formed near the antiferromag- 
netic quantum critical point (QCP)^ as the pressure or in 
other cases the concentration of a component is varied. 
However, this theory seems to be insufficient to explain 
the temperature-pressure phase diagram of CeCu2Ge2 or 
CeCu2Si2, where a superconducting dome with enhanced 
transition temperature is located far away from the an- 
tiferromagnetic critical point. Since this discovery, these 
compounds have drawn much attention both experimen- 
tally and theoretically. It has been argued that this phe- 
nomenon is related to the critical valence fluctuations of 
the Ce ionsr'"— that is, to the existence of a second QCP. 

The simplest model of heavy-fermion compounds is the 
periodic Anderson model (PAM) It is known, however, 
that the mixed-valence regime appears always in this 
model as a smooth crossover, and valence fluctuations do 
not become critical for any choice of the parameters. A 
local Coulomb interaction between the conduction and 
localized electrons is needed for the appearence of a 
sharp transition and critical valence fluctuations i^*^ The 
Hamiltonian of this extended periodic Anderson model 
can be written in the following form: 

(1) 



where the notations are standard. After Onishi and 
Miyake's pioneering work)^ this model has been inves- 
tigated by several modern techniques recently, including 
density matrix renormalization group,— dynamical mean 
fieldr^"— variational calculations f^ii^ projector based 
renormalization approach^, and fluctuation exchange 
approximation<2i It has been found that a first-order va- 
lence transition and a QCP may appear due to Udf. It is 
worth noting that a first-order transition from Mott in- 
sulator to Kondo insulator has been found^^ in a model 
with a more general Hamiltonian, too, including Hund's 
coupling and interaction between c?-electrons. The pre- 
vious calculations focused mainly on the infinite or large 
U f cases. Little is known on the properties of the sys- 
tem near the critical point, or how the QCP depends 
on the parameters of the model. Our main goal in this 
paper is to further examine this model with mean-field 
theory, Gutzwiller variational method, and exact diag- 
onalization. We will study the properties of the model 
defined in ([1]) for a wide range of parameters in the half- 
filled case. 



The setup of the paper is as follows. In Sec. II, we 
perform a mean-field calculation to demonstrate in the 
simplest way how Udf affects the intermediate valence 
regime. In Sec. Ill, the variational approach is intro- 
duced, which is based on the Gutzwiller wave function. 
We analyze the quantum critical behavior and the disap- 
pearance of the Kondo regime at a triple point. More- 
over, we construct the £/ — Udf phase diagram. In Sec. 
IV, we carry out exact diagonalization to investigate the 
model in one dimension and compare the results with 
that of mean-field theory and the Gutzwiller approach. 
Finally, in Sec. V, our conclusions are presented. 
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II. MEAN-FIELD CALCULATION 

First of all we will study the problem by mean-field 
methods. In this approach some kind of order has to be 
assumed. We performed calculations by assuming sev- 
eral possibilities: a) the system is paramagnetic; b) it 
possesses magnetic order, which can be ferromagnetic or 
of spin-density-wave type. The state with lowest energy 
is accepted as the ground state. According to our calcu- 
lations, the mean-field equations always have a paramag- 
netic solution, but in the Kondo regime, where localized 
moments are present, ordering of the moments leads to 
the lowering of the energy. Since the ferromagnctically 
ordered state has higher energy than the spin-density- 
wave state, the latter will be the subject of our investiga- 
tions. We assume a simple cubic lattice in the following 
calculations, which can be partitioned into two sublat- 
tices {A and B). We expect that in the broken symmetry 
phase the electrons are ordered on the two sublattices in 
an alternating fashion, that is: 



•J 



Jq,rRj 



(2) 



where is the magnetization of the sublattice and <7q = 
7r/a(l,l,l), so that QgRj = l-nn on sublattice A and 
q^Rj — i^n + l)7r on sublattice B {n is an integer), 
and nj is the average number of /-electrons per site. 
The values of n/ and m/ need to be determined self- 
consistcntly. Similar oscillation can be assumed for the 
d-electrons, 



nd + crmde 



(3) 



although nid will not appear explicitly in the calculations. 
The mean-field Hamiltonian is 



NUf , 2 2\ 

— — [jij- mj) - NUdfUfUd, 



(4) 



where the fe sum extends over the whole Brillouin zone 
of the simple cubic lattice. Due to the assumed mag- 
netic ordering, the size of the Brillouin zone is reduced 
to half of its original size. In order to restrict the k sum 
to the magnetic Brillouin zone, we split the original sum 
into two parts by introducing the operators dl.^^^^ and 

fk+qory- suppose that the dispersion relation pos- 
sesses the nesting property: 



+ Qo) = -£d(fc), 



(5) 



which is valid in a tight-binding model with nearest 
neighbour hopping. This fixes the zero of the energy 
scale. Then the mean-field Hamiltonian can be rewritten 
in Bloch representation: 



n 



AF 



fka 
V fk+qoa 



t 



n{k,a) 



dk+qocr 
Ik a 
V fk+qo<T / 
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mj) - NUdfUfUd, 



(6) 



where the prime denotes that the summation is carried 
out over the magnetic Brillouin zone, and 



( Uk) 

Uk) 
V 
\ V 



V ^ \ 

V 



(7) 



where ^/ = £/ -f Ufnf/2 + UdfUd, U{k) = ed{k) + UdfUf 
and £,d{k) — ~£d{k) + Udfnf. This Hamiltonian can be 
diagonalizcd by the unitary transformation T{ka), which 
is a real matrix in our case, leading to 

^AF = E ' [Eaik)AUk. + Et{k)BlBka 

k.a 

+Ecik)Ci^Cka + Edik)Dl^Dka 
_^[EL ("-/ ^ "^z) ^ NUdfUfUd, 



(8) 



where 



( 



A^ rT /^T rit 
^fccr ^ka '-ka ^ka 



^fccr ^k+q„a fka -f I 



t 



T\ka). 



(9) 



The diagonalization is done numerically for each fe 
value and wc sort the eigenvalues in increasing order 
[Ea{k) < Eb{k) < Ec{k) < Ed{k)]. Moreover, the eigen- 
values have to be determined iteratively, since n/, to/ 
and Ud appearing in 'H{k,<7) have to satisfy a self- 
consistency condition. This condition can easily be for- 
mulated in the half-filled case, where — as it will be dis- 
cussed below — the two lower bands with dispersion Ea{k) 
and Eb{k) are fully occupied and the two higher lying 
bands are empty. The conditions of self-consistency for 
Uf and nif are 



^ E ' [^31 {k<T) + Ti^ (ka) + Tl (ka) +Tlika] 

k,C7 



(10) 



k 



(11) 
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and the total energy is 

k,(T 

-^{n}-m))-NUafnfna. (12) 

Before evaluating the self-consistency equations, we re- 
turn to the problem of the eigenvalue equations. Analytic 
expressions can be given at the symmetric point, and the 
obtained four bands for Udf = become simply 



± 



(elik)-Ujm}) +W^[el{k) + Ujm}) 



(13) 



The band structure is displayed in Fig. [T] for the one- 
dimensional tight-binding case for a special, but non- 
symmetric choice of the parameters. In higher dimcn- 




FIG. 1: The dispersion curves of the states diagonalizing 
the mean-field Hamiltonian for Uf/W = 0.5, Udf /W = 0, 
Ef/W = -0.1, V/W = 0.1. 

sions, the system docs not necessary remain insulating if 
the gap opens at different energies in different points of 
the zone boundary. 

The self-consistency equations always have a param- 
agnetic solution, ruf = 0, and in certain cases they have 
a magnetic solution, |m/| > 0. We compare the energies 
of both solutions and accept that one. which has lower 
energy. It is worth noting that the polarization of the d- 
electrons, m^, is also nonzero and its sign is opposite to 
TO f in the magnetic solution (though its value is smaller 
than TO/ by an order of magnitude) . 

In the actual calculations we do not work in fc-spacc. 
The summations over k are carried out by assuming a 
constant density of states p{e) = 1/W in the interval 
e e [-W/2,W/2]. The same p(e) wiU be used in the 
variational calculation, too. The numerical results for 
Uf are shown in Fig. [2] for different values of Udf- As 




Ef/W 

FIG. 2: The /-level occupancy as a function of e/ obtained by 
mean-field theory. The solid, dashed and dotted lines belong 
to Udf/W = 0;0.4; 1, respectively, V/W = 0.1, Uf/W = 3 in 
all cases. 



long as the /-level is nearly fully occupied or empty, 
the paramagnetic solution is favorable, while when the 
occupancy is nearly 1, the magnetic solution has lower 
energy. For Udf/W = 0, the /-level occupancy is contin- 
uous, although there is a discontinuity in its derivative 
at the point where the paramagnetic solution switches 
to magnetic or vice versa. This happens in Fig. [2] at 
Ef/W sa —3.45 and 0.45. For increasing Udf values, 
first the mixed-valence regime narrows, then a jump, a 
first-order transition appears in the /-level occupancy be- 
tween the paramagnetic and the magnetic solution and 
the Kondo regime shrinks rapidly. This first-order transi- 
tion has already been found by other calculations J^ii^ As 
is seen here, this simple mean-field approach can also ac- 
count for it. The mean-field theory gives a critical value 
of U^j:/W « 0.26 for V/W = 0.1, Uf/W = 3. 

At a certain value of Udf above U^^- the magnetic solu- 
tion is stable in a single point, £/ = —U //2. In this point 
the paramagnetic solutions with « 2 and ny « 0, and 
the magnetic solution with n/ = 1 have the same energy, 
that is, three states coexist. This is a triple point, since 
three first-order transition lines meet here. For larger 
Udf values, beyond the triple point, a so-called valence 
skipping occurs, that is, the valence state « 1 is miss- 
ing, since a direct first-order transition takes place from 
71/ « 2 to n/ ~ 0. It is interesting to note that so far 
the valence skipping (which was observed in several com- 
punds, for example BaBiOs) has been attributed to the 
presence of a negative Uf^ As it is demonstrated here, 
large enough Udf can also lead to valence skipping, even 
if Uf > 0. The mean-field theory gives U^f^^^yW w 1.75 
for V/W = 0.1, Uf/W = 3. Note that this is in good 
agreement with the result U^f^^^° ~ C///2 + VF/4 that will 
be obtained by the Gutzwiller method. 

This picture remains valid even when Uf is small, but 
finite, compared to the bandwidth. Although no plateau 
with w 1 is formed, that is, there is no Kondo-rcgime, 
a stable magnetic solution is found near the symmetric 
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point, and this regime is bounded by first-order transi- 
tions above a finite U^j- Note that for Uf/W = there 
is no triple point due to the lack of magnetic solution. 
However, a direct first-order transition appears between 
the paramagnetic states with solution between nj k 2 
and n/ « above a certain value of Udf- 

Although this theory gives a qualitatively good de- 
scription of the first-order transition, there are several 
problems with it. Firstly, we had to assume a magnetic 
order of spin-density-wave type, while no long-range or- 
der is expected in the Kondo regime. Second, the valence 
susceptibility, which is defined by 



relations between d- and /-electrons, for example: 



XV 



def ' 



(14) 



is not a continuous function, as is mentioned above, even 
below U^jr/W, where n/ is continuous. The mean-field 
approach does not provide us with a critical point, where 
Xv would diverge. In order to find a QCP and to inves- 
tigate its properties, we need a more accurate calcula- 
tion. This is done in the next section by the variational 
method. 



III. VARIATIONAL CALCULATION 

In what follows we generalize the variational approach 
used in [Ref. [Hll and summarize briefly the main steps 
of the calculation. We restrict ourselves to the param- 
agnetic case, that is, the number of up-spin and down- 
spin electrons are assumed to be equal locally, too. As it 
will be pointed out, the quantum criticality and the first- 
order transition appear without any further assumptions, 
in contrast to the mean-field calculation. The trial state 
is expressed in terms of Gutzwiller projectors: 



!*) = p{f'd')p{fd')p{fd^)p{fd')p{fd^)\^o), 



(15) 



where 



i*o) = n n + J io) (16) 



contains the mixing amplitudes u{k) and v{k) as vari- 
ational parameters, and the sum over k extends over 
the whole Brillouin zone. The Gutzwiller projectors 
P{f"d^), where a and /3 denote the /- and d-electron 
numbers respectively, act on the on-site electron config- 
urations defined by their arguments, making that config- 
uration less probable. For example: 



Pif'd 



4Ai 



(17) 



is the standard Gutzwiller-projector for two /-electrons 
on the same site. The other ones take into account cor- 



[^-Vif'd')) 



g,)\KA}^ (18) 
p{m = n {i - (1 - '?(/''^')) KAiKA} ■ 

g 

(19) 

The remaining projectors arc defined straightforwardly. 
Besides the mixing amplitudes, we have five variational 
parameters, v{pd^), viPd^), viPd^), viPd'), v(.Pd°), 
controlled by Udf and U / . The tedious procedure of op- 
timization is omitted here. Performing the optimiza- 
tion with respect to the mixing amplitudes using the 
Gutzwiller approximation we obtain 



£ = 



-E 



qdSdik) + £f - \l [qd£d{k) - £/] ^ + W"^ 



+(£/ - ef)nf + Ufi^ifd^) 

+Udf [Mfd') + 2 {u{fd') + v{pd^)) + J^ipd^)] 

(20) 

for the ground-state energy density, where v{pd^) is the 
density of doubly occupied /-sites. The other iy{f"d^) 
quantities denote the corresponding densities of the f"d^ 
configurations, e.g.: 



'at 



1 - h 



gl 



^gi 



^gt 



(21) 



V = Vy^/qjqd is the renormalized hybridization, qf and 
qd are the kinetic energy renormalization factors for the 
/- and d-electrons, respectively. Their analytic forms are 
now much longer than in our previous paper )^ therefore 
they are not presented here. However, it is remarkable 
that the renormalized hybridization can still be written 
as the square root of g/ and qd in the presence of Udf, 
too. Furthermore, if is the quasiparticle energy level of 
/-electrons; it has the same form as in [Ref. [25|. It pro- 
vides a self-consistency equation for nf. The summation 
over k in Eq. (PH)) is carried out with a constant density 
of states, p{e) = 1/W, in the interval e € [-W/2, W/2]. 
During the optimization process, the r]{f"d^) variational 
parameters are expressed as functions of the quantities 
i'{f°'d^), therefore the actual optimization can be done 
with respect to these parameters. All in all, the energy 
density given in (|20p has to be optimized for Uf and 
i^(/2dO), i^iPd'), iy{pd^), u{pd^), v{pd^). The evalu- 
ation of these equations could be done only numerically. 

We first address what happens in the mixed-valence 
regime. As Udf is switched on, the mixed-valence regimes 
tend to be sharper and sharper. This can be character- 
ized by the valence susceptibility defined in Eq. (|14l) and 
displayed in Fig. [3] for different values of Udf - Note that 
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in the half-filled case, this function is symmetric to the 
point £/ = —Uf/2, therefore it is sufficient to investi- 
gate the critical behavior in the regime < riy < 1. 
In what follows, we focus on this regime, if not men- 
tioned otherwise. It is found, in agreement with other 




FIG. 3: The valence susceptibility as a function of the /-level 
energy for (7//VK = 3 and V /W = Q.l. The solid, dashed and 
dotted lines correspond to Udf/W — 0; 0.3 and 0.4, respec- 
tively. 

calculations r^ii^ that xv diverges for a certain value of 
U^f and two values of e^- related by the symmetry with 
respect to —Uf/2. These two points are identified as the 
QCPs. Following the maximum values of xvy the posi- 
tion of the QCP can be determined. We found that xv 
diverges as 



1 



Xv\u^j=U2f 



1 



(22) 
(23) 



This power-law behavior is valid for every choice of the 
parameters we used in our calculations, indicating uni- 
versality. It is interesting to investigate how the position 
of the QCP depends on Uf and V. In Fig. |3]the critical 
U^f and (for the QCP in the < < 1 regime) are 
shown as a function of Uf for a fixed V. We found that 
(i) even for C// = there exists a critical point; (ii) U^f 
and e'j vary monotonically as Uf increases; (iii) both of 
them saturate as [// reaches the value above which there 
exists a stable Kondo regime (see Fig. 3. in [Ref.[25|). On 
the contrary their dependence on V is remarkable. These 
values are shown in Fig. [S] Firstly, we mention, that for 
V/W 0, U^f tends to a nonzero value {U^f/W « 0.17), 
which indicates that there has to be a finite value of Udf 
even for weak hybridizations to obtain a valence transi- 
tion. The critical position of the /-level decreases linearly 
with increasing hybridization. The critical value of the 
occupancy of the /-level increases from = at ^ = 
and saturates as soon as reaches the bottom of the 
conduction band. Roughly at the same mixing V, the 
slope of the U^j — V curve shows a substantial change. 
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FIG. 4: Panel (a) shows the critical value of Udf where the 
QCP appears as a function oiUf. Panel (b) shows the critical 
value of £/, V/W = 0.1 in all cases. 



For larger values of Udf, two subsequent first-order 
transitions — from Uf 2 to nf I (Kondo regime) and 
from Ti^wltorij^sa — take place as is varied. Their 
positions are symmetric with respect to —Uf/2. This is 
confirmed by the fact that near the transition a hysteresis 
is observed, that is, there is a narrow range of £/, where 
two solutions of the variational equations coexist. There- 
fore the transition line is identified from the ground-state 
energy, where the energies of the different configurations 
are equal. 

For even larger values of Udf, the width of Kondo 
regime decreases and at J7^"^'° it ends in a triple point. 
At the triple point, the energy of the Kondo- like state 
becomes equal to the energy of the states with w 2 
and Uf ~ 0; therefore here three different states co- 
exist. We found that the triple point is located at 



triple 



-Uf/2 and U^ 



triple 
df 



Uf/2 + W/A, if there 
is a Kondo plateau. For small Uf (includig Uf = 0), 
when there is no Kondo plateau, the numerical results 
can be fitted to C/]}'p'° w C///2 + W/3. In [Ref. [11] using 
DMFT it was found that the Kondo regime is stable for 
Udf ^ Uf/2. Our result is in agreement with this. The 
mean-field theory gives similar results except for Uf = 0, 
where the triple point does not exist. Now we can draw 
the Sf — Udf phase diagram. The results are shown in 
Fig. [51 using a color code, for two different values of 
the hybridization and demonstrates our statements de- 
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FIG. 5: (Color online) Panel (a) shows the critical value of Ud] 
(•) where the QCP appears as a function of the hybridization. 
The red dotted and blue dashed lines are linear fits to the be- 
ginning and to the end of the data respectively. Panel (b) and 
(c) show the critical values of £/ and critical n/ respectively, 
Uj jW = 3 in all cases. 



scribed above. The figure demonstrates that the interval 
of £/, where first-order transition occurs, is shortened for 
larger values of the hybridization. 



IV. COMPARISON WITH EXACT 
DIAGONALIZATION 

As a check of the reliability of the mean-field and vari- 
ational results, which are known to be accurate in infi- 
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FIG. 6: (Color online) The upper panel shows the £/ — Ud; 
phase diagram for UjjW — 3 and V/W — 0.1, while the 
lower one for Uf /W = 3 and V/W = 0.2. The solid white 
lines denote the first-order transition lines. The QCPs are 
marked by red squares. The point where the white lines meet 
is the triple point (see the text). 



nite dimensions, we performed exact diagonalization on a 
one-dimensional chain. Due to the limitation to relatively 
short chains containing six sites, wc do not expect to find 
critical behavior in this calculation. However, some other 
features of the effect of Udf might be observable. The 
comparison is shown in Fig. [7] for a special choice of the 
parameters. The width of the Kondo plateau is the same 
using all the three methods, therefore its shrinking due 
to Udf is not an artifact of the Gutzwiller approximation 
or the mean-field treatment. Furthermore, by increasing 
Udf, the intermediate valence regime becomes narrower 
in the exact diagonalization, too, although there is nat- 
urally no sharp valence transition. 

The mean-field theory and the Gutzwiller approach 
yield surprisingly close results for Udf/W = and 
1. However, the mean-field results show a jump in 
the /-level occupancy at such small values of Udf, 
where the Gutzwiller method still gives a continuous 
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FIG. 7: The /-level occupancy as a function of e/. The black 
curves are obtained from the Gutzwiller method, while the 
gray ones are calculated from exact diagonalization and the 
symbols are the results of the mean-field calculation. The 
solid, dashed and dotted lines (and the symbols •, A, ■) 
belong to Udf /W = 0; 0.4; 1 respectively, V/W = 0.1 and 
Uf/W = 3 in all cases. 

change oi Uf. The estimated eritieal value (U^j^) from 
the mean-field theory is signifieantly smaller than that 
from the Gutzwiller method. It is worth emphasizing 
that the jump in the mean-field results originates from 
the paramagnetic- magnetie transition. In contrast, the 
Gutzwiller method gives a valence transition between 
paramagnetic states. Both methods result in a triple 
point for a certain value of Udf, and we found that both 
of them gives C/^^'P''' KUf/2 + W/4. For small values of 
Uf the scenario is the same in both methods as in the 
large Uf case, however, the values of U^"^^° and U^^ are 
somewhat different. The only exception is [// = 0, where 
there is no triple point in the mean-field calculation due 
to the missing of a stable magnetic solution. Here a di- 



rect first-order transition takes place between the nearly 
fully occupied and nearly empty /-levels. 



V. CONCLUSIONS 

We have performed a mean-field calculation, a varia- 
tional calculation using the Gutzwiller method, and exact 
diagonalization for the extended PAM, where an addi- 
tional local Coulomb interaction between the d- and /- 
electrons has been included. Both the mean-field theory 
and the Gutzwiller method have resulted in two subse- 
quent first-order valence transitions as the position of the 
/-level is varied above a critical value of Udf- Thus two 
QCPs appear in the £/ — Udf plane. There is only one 
QCP in the mean-field calculation at Uf/W = 0, while 
the scenario does not change in the Gutzwiller-method. 
We have analyzed variationally the critical behavior as a 
function of hybridization, the bare /-level energy, and U f , 
and have drawn the £/ — Udf phase diagram. It has been 
pointed out that the Kondo regime shrinks by increasing 
Udf, and ends in a triple point. For even larger values of 
Udf a direct first-order valence transition takes place from 
rij sa 2 to njsaO. This can be interpreted as valence 
skipping, which so far has been attributed to the pres- 
ence of a negative Uf. We find it for Uf > 0, when Udf is 
large enough. The shrinking of the Kondo regime and the 
narrowing of the intermediate valence regime have been 
confirmed by exact diagonalization, although naturally, 
no sharp valence transition is found in finite chains. 
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